***File used to compute results and generate tables for Murphy "Sobriety, Social Capital, and Village Network Structures"



cd "User File Location"

clear all
macro drop _all

*CUSTOM DO FILES

///User must manually install nwcommands-ado from http://www.nwcommands.org or by using the search command and selecting "install".

/*
ssc install balancetable
ssc install estout
*/

use "alcdatacleanednets", clear
sort hhid indvid

drop if attrited==1
drop if purecontrol==1

gen treat1_2=.
replace treat1_2=1 if treatment==1
replace treat1_2=0 if control==1
gen treat1_3=.
replace treat1_3=1 if treatment==1
gen treat2_3=.
replace treat2_3=1 if control==1

label var blalcoholnumber "Baseline: Alcoholic drinks past week$ \S $"
global balance iteso luhya blage tlu blalcoholnumber blagworkownfarmhours blagworkdifffarmhours blnonagwageshours blhouseworkhours mathproblem numnearbypeers  avgdistdemean landareaacres assetindx nightsaway hhsize readenglish readswahili distvillagemediankm farmer yrseduc multiplespouses yearslivedinvillage
global balancesplitb landareaacres assetindx nightsaway hhsize readenglish readswahili distvillagemediankm farmer yrseduc multiplespouses yearslivedinvillage
global balancesplita iteso luhya blage tlu blalcoholnumber blagworkownfarmhours blagworkdifffarmhours blnonagwageshours blhouseworkhours  distvillagemediankm avgdistdemean numnearbypeers 
global balancesplitb landareaacres assetindx nightsaway hhsize readenglish readswahili mathproblem farmer yrseduc multiplespouses yearslivedinvillage



///Balance tables
hotelling $balance, by(treat1_2)
local F=((r(N)- r(k)-1)/((r(N)-2)*r(k)))*r(T2)
local fpvalue1=round(Ftail(r(k), r(df), `F'),0.001)

balancetable (mean if blgroup!=3) (mean if blgroup==1) (mean if blgroup==2) (diff treatment if blgroup!=3) $balancesplita using "groupbalance2nopca.tex", varlabels replace ctitles("Overall Mean" "Mean Treat" "Mean Control" "Treat $ - $ Control")
balancetable (mean if blgroup!=3) (mean if blgroup==1) (mean if blgroup==2) (diff treatment if blgroup!=3) $balancesplitb using "groupbalance2nopcb.tex", varlabels replace ctitles("Overall Mean" "Mean Treat" "Mean Control" "Treat $ - $ Control") prefoot("Hotelling test (p-value)& & & & `fpvalue1' \\") 



use "dyadicdatasetforestimation", clear
global villagenumber villagenum2 villagenum3 villagenum4 villagenum5 villagenum6 villagenum7 villagenum8 villagenum9 villagenum10 villagenum11 villagenum12 villagenum13 villagenum14 villagenum15 villagenum16 villagenum17 villagenum18 villagenum19 villagenum20 villagenum21 villagenum22 villagenum23 villagenum24 villagenum25 villagenum26 villagenum27 villagenum28 villagenum29 villagenum30
global covariateslimited iteso luhya tlu mathproblem readenglish readswahili multiplespouses hhsize blage yrseduc assetindx nightsaway blalcoholnumber distvillagemediankm landareaacres farmer
global covariateslimitedmatch itesomch luhyamch tlumch mathproblemmch readenglishmch readswahilimch multiplespousesmch hhsizemch blagemch yrseducmch assetindxmch nightsawaymch blalcoholnumbermch distvillagemediankmmch landareaacresmch farmermch

global sumvars1
global diffvars1
local num1 : word count $covariateslimited
rename infopublichealth health

forvalues x=1/`num1'{
local z : word `x' of $covariateslimited
local y : word `x' of $covariateslimitedmatch
gen `z'sum=`z'+`y'
gen `z'diff=`z'-`y'
global sumvars1 $sumvars1 `z'sum
global diffvars1 $diffvars1 `z'diff 
}


gen distancekm3=distancekm^3
gen distancekm4=distancekm^4


///Link Predictions
logit knowperson treattreat treatcontrol controltreat stddist malemale malefemale femalemale $sumvars1 $diffvars1 i.villagenumber enumeratormale1_1 enumeratormale1_2 enumeratormale2_1 enumeratorfemale1_1 enumeratorfemale1_2 enumeratorfemale2_1 if indvid!=indvidmch, cluster(indvid)


predict pred_know if indvid!=indvidmch

gen bipred_know=1 if pred_know>.64&indvid!=indvidmch
replace bipred_know=0 if bipred_know==.&indvid!=indvidmch
gen correctpred=1 if bipred_know==1&knowperson==1
replace correctpred=1 if bipred_know==0&knowperson==0
replace correctpred=0 if bipred_know==1&knowperson==0
replace correctpred=0 if bipred_know==0&knowperson==1
tab correctpred






logit advice treattreat treatcontrol controltreat stddist spouse malemale malefemale femalemale $sumvars1 $diffvars1 i.villagenumber enumeratormale1_1 enumeratormale1_2 enumeratormale2_1 enumeratorfemale1_1 enumeratorfemale1_2 enumeratorfemale2_1  if indvid!=indvidmch, cluster(indvid)
predict pred_advice if indvid!=indvidmch

eststo log1: margins, dydx(treattreat treatcontrol controltreat stddist malemale malefemale femalemale spouse) predict(pr) post
estadd local type "Logit"
estadd local regressors "Yes"
estadd local fe "Yes"

gen bipred_advice=1 if pred_advice>.265&indvid!=indvidmch
replace bipred_advice=0 if bipred_advice==.&indvid!=indvidmch
gen correctpredadvice=1 if bipred_advice==1&advice==1
replace correctpredadvice=1 if bipred_advice==0&advice==0
replace correctpredadvice=0 if bipred_advice==1&advice==0
replace correctpredadvice=0 if bipred_advice==0&advice==1
tab correctpredadvice 


logit trust treattreat treatcontrol controltreat stddist spouse malemale malefemale femalemale $sumvars1 $diffvars1 i.villagenumber enumeratormale1_1 enumeratormale1_2 enumeratormale2_1 enumeratorfemale1_1 enumeratorfemale1_2 enumeratorfemale2_1  if indvid!=indvidmch, cluster(indvid)
predict pred_trust if indvid!=indvidmch

eststo log2: margins, dydx(treattreat treatcontrol controltreat stddist malemale malefemale femalemale spouse) predict(pr) post
estadd local type "Logit"
estadd local regressors "Yes"
estadd local fe "Yes"


gen bipred_trust=1 if pred_trust>.215&indvid!=indvidmch
replace bipred_trust=0 if bipred_trust==.&indvid!=indvidmch
gen correctpredtrust=1 if bipred_trust==1&trust==1
replace correctpredtrust=1 if bipred_trust==0&trust==0
replace correctpredtrust=0 if bipred_trust==1&trust==0
replace correctpredtrust=0 if bipred_trust==0&trust==1
tab correctpredtrust


logit speakdaily treattreat treatcontrol controltreat stddist spouse malemale malefemale femalemale $sumvars1 $diffvars1 i.villagenumber enumeratormale1_1 enumeratormale1_2 enumeratormale2_1 enumeratorfemale1_1 enumeratorfemale1_2 enumeratorfemale2_1  if indvid!=indvidmch, cluster(indvid)
predict pred_daily if indvid!=indvidmch
eststo log3: margins, dydx(treattreat treatcontrol controltreat stddist malemale malefemale femalemale spouse) predict(pr) post
estadd local type "Logit"
estadd local regressors "Yes"
estadd local fe "Yes"



gen bipred_daily=1 if pred_daily>.30&indvid!=indvidmch
replace bipred_daily=0 if bipred_daily==.&indvid!=indvidmch
gen correctpreddaily=1 if bipred_daily==1&speakdaily==1
replace correctpreddaily=1 if bipred_daily==0&speakdaily==0
replace correctpreddaily=0 if bipred_daily==1&speakdaily==0
replace correctpreddaily=0 if bipred_daily==0&speakdaily==1
tab correctpreddaily


logit workedwith treattreat treatcontrol controltreat stddist spouse malemale malefemale femalemale $sumvars1 $diffvars1 i.villagenumber enumeratormale1_1 enumeratormale1_2 enumeratormale2_1 enumeratorfemale1_1 enumeratorfemale1_2 enumeratorfemale2_1  if indvid!=indvidmch, cluster(indvid)
predict pred_workedwith if indvid!=indvidmch

eststo log4: margins, dydx(treattreat treatcontrol controltreat stddist malemale malefemale femalemale spouse) predict(pr) post
estadd local type "Logit"
estadd local regressors "Yes"
estadd local fe "Yes"

gen bipred_workedwith=1 if pred_workedwith>.255&indvid!=indvidmch
replace bipred_workedwith=0 if bipred_workedwith==.&indvid!=indvidmch
gen correctpredworkedwith=1 if bipred_workedwith==1&workedwith==1
replace correctpredworkedwith=1 if bipred_workedwith==0&workedwith==0
replace correctpredworkedwith=0 if bipred_workedwith==1&workedwith==0
replace correctpredworkedwith=0 if bipred_workedwith==0&workedwith==1
tab correctpredworkedwith



logit health treattreat treatcontrol controltreat stddist spouse malemale malefemale femalemale $sumvars1 $diffvars1 i.villagenumber enumeratormale1_1 enumeratormale1_2 enumeratormale2_1 enumeratorfemale1_1 enumeratorfemale1_2 enumeratorfemale2_1  if indvid!=indvidmch, cluster(indvid)
predict pred_health if indvid!=indvidmch

eststo log5: margins, dydx(treattreat treatcontrol controltreat stddist malemale malefemale femalemale spouse) predict(pr) post
estadd local type "Logit"
estadd local regressors "Yes"
estadd local fe "Yes"

gen bipred_health=1 if pred_health>.368&indvid!=indvidmch
replace bipred_health=0 if bipred_health==.&indvid!=indvidmch
gen correctpredhealth=1 if bipred_health==1&health==1
replace correctpredhealth=1 if bipred_health==0&health==0
replace correctpredhealth=0 if bipred_health==1&health==0
replace correctpredhealth=0 if bipred_health==0&health==1
tab correctpredhealth



logit closefriend treattreat treatcontrol controltreat stddist spouse malemale malefemale femalemale $sumvars1 $diffvars1 i.villagenumber enumeratormale1_1 enumeratormale1_2 enumeratormale2_1 enumeratorfemale1_1 enumeratorfemale1_2 enumeratorfemale2_1  if indvid!=indvidmch, cluster(indvid)
predict pred_closefriend if indvid!=indvidmch

eststo log6: margins, dydx(treattreat treatcontrol controltreat stddist malemale malefemale femalemale spouse) predict(pr) post
estadd local type "Logit"
estadd local regressors "Yes"
estadd local fe "Yes"

gen bipred_closefriend=1 if pred_closefriend>.225&indvid!=indvidmch
replace bipred_closefriend=0 if bipred_closefriend==.&indvid!=indvidmch
gen correctpredclosefriend=1 if bipred_closefriend==1&closefriend==1
replace correctpredclosefriend=1 if bipred_closefriend==0&closefriend==0
replace correctpredclosefriend=0 if bipred_closefriend==1&closefriend==0
replace correctpredclosefriend=0 if bipred_closefriend==0&closefriend==1
tab correctpredclosefriend


esttab log1 log2 log3 log4 log5 log6  using "logitforpred.tex", label replace obslast star(* 0.10 ** 0.05 *** 0.01) se nomtitles mgroups("Advice from $ j $" "Trust $ j $" "Speak with $ j $ daily" "Worked with $ j $" "Health Info from $ j $" "Close Friend with $ j $"  , pattern(1 1 1 1 1 1) prefix(\multicolumn{@span}{c}{) suffix(}) span erepeat(\cmidrule(lr){@span})) substitute(\_ _) ///
scalars("fe Village/Enumerator FEs" "regressors Additional symmetric regressors$ \ddagger $") b(2) sfmt(2) eqlabels(none) nonotes nogaps noomitted order(treattreat treatcontrol controltreat stddist spouse malemale malefemale femalemale) drop()







save "dyadicdatanetworkanalysis.dta", replace






///Centrality Analysis///
nwclear
use "dyadicdatanetworkanalysis.dta", clear
sort villagenumber indvid indvidmch

global bipred advice trust daily workedwith health closefriend
local r=1

nwclear



forvalues e=1/6{
local asx: word `e' of $bipred

foreach y in 1 2 3 5 7 8 10 11 12 14 16 17 18 19 21 24 26 27 28 30{
nwclear
use "dyadicdatanetworkanalysis.dta", clear
sort villagenumber indvid indvidmch
keep if villagenumber==`y'
drop if attrited==.
nwset indvid indvidmch bipred_`asx', name(network_`asx'`y') edgelist keeporiginal
rename _nodeoriginal indvid

merge 1:1 indvid using "alcdatacleanednets"
drop if _merge!=3
drop _merge

nwsave network`asx'`y', replace


nwclear
use "dyadicdatanetworkanalysis.dta", clear
sort villagenumber indvid indvidmch

gen n=_n
keep if villagenumber==`y'
drop if attrited==.
gen samehh=1 if hhid==hhidmch
replace samehh=0 if hhid!=hhidmch
replace samehh=0 if indvid==indvidmch
nwset indvid  indvidmch samehh, name(network_samehh`y') edgelist keeporiginal
rename _nodeoriginal indvid
nwsave networkhh`y', replace

}

clear
nwclear


global villagename Akibui Kongurakol Akiriamet Apokor Matisi Bokokholo Khasolo Abulei Emusioma Misikhu Malaha Sango Sirisia Naika Lunakwe Namisi Shipalo Shiraha Shirali Emuchelela Imanga Emokokho Milimani Kokare Kiriko Kocholia Lwanda EmukokhoA KapijanB Kakoit
global graphtitle Receive_Ag_Advice Trust Speak_Daily Worked_with Receive_Health_Info Close_Friend 
foreach y in 1 2 3 5 7 8 10 11 12 14 16 17 18 19 21 24 26 27 28 30{
clear 
nwclear
nwuse networkhh`y'
nwuse network`asx'`y'

local z: word `y' of $villagename
local title: word `e' of $graphtitle

gen n=_n
egen maxn=max(n)
replace maxn=maxn-1

nwcloseness network_`asx'`y', unconnected(max)

nwdegree network_`asx'`y'

nwbetween network_`asx'`y'

nwdegree network_`asx'`y', valued
sort indvid

/// Network figures (note - nwplot does not necessarily create identical graph each time code is run, but graphs should appear similar) 
nwplot network_`asx'`y', edgecolor(network_samehh`y', edgecolorpalette(gray gold)) title(`z' Village Network: `title') label(indvid) arcstyle(straight) legend(label(1 "Male") label(2 "Female") label(3 "Non-HH link") label(4 "HH link")) labelopt(mlabsize(tiny)) iterations(100) color(blfemale, colorpalette(blue red)) layout(mds) nodefactor(.5) edgefactor(.5) arrowfactor(.1) arrowbarbfactor(.1)
graph export networkvillage`y'graph`r'.pdf, replace


global outcomevars _in_degree _out_degree _close _between 
foreach var of varlist $outcomevars{
gen `var'std=.
sum `var' , detail 
replace `var'std=((`var'-`r(mean)')/`r(sd)') 
}
save Stub`y', replace
}

nwclear
use Stub1, clear

foreach y in 2 3 5 7 8 10 11 12 14 16 17 18 19 21 24 26 27 28 30{
append using Stub`y'

}
keep indvid n-_betweenstd
save CentralityMeasures, replace
foreach xe of varlist indvid-_betweenstd{
rename `xe' `xe'mch
}
save CentralityMeasuresmch, replace





///Centrality regressions///

use "alcdatacleanednets", clear
drop if indvid==.
drop if purecontrol==1
merge 1:1 indvid using CentralityMeasures
drop _merge
gen treatmentmale=treatment*blmale
gen controlmale=control*blmale
global villagecontrols blmale blage blalcoholnumber landareaacres assetindx nightsaway hhsize readenglish distvillagemediankm farmer yrseduc multiplespouses enumeratormale1 enumeratorfemale1 i.villagenumber
label var treatmentmale "Treatment assignment $ \times $ male"
label var distvillagemediankm "Distance to village center (kilometers)"

local q=1
global outcomevars1 _in_degreestd _out_degreestd _closenessstd _betweenstd 


foreach t of varlist $outcomevars1{
eststo reg`q'_`r': reg `t' treatment blmale enumeratormale1 enumeratorfemale1 i.villagenumber, vce(bootstrap, reps(1000) seed(2345) cluster(hhid))
estadd local type "ITT"
estadd local fixed "Yes"
estadd local regressors "No"


eststo reg`q'a_`r': reg `t' treatment $villagecontrols, vce(bootstrap, reps(1000) seed(2345) cluster(hhid))
estadd local type "ITT"
estadd local regressors "Yes"
estadd local fixed "Yes"


eststo reg`q'b_`r': reg `t' treatment treatmentmale $villagecontrols, vce(bootstrap, reps(1000) seed(2345) cluster(hhid))
estadd local type "ITT"
estadd local regressors "Yes"
estadd local fixed "Yes"


local ++q

}

esttab reg1_`r' reg1a_`r' reg1b_`r' reg2_`r' reg2a_`r' reg2b_`r' reg3_`r' reg3a_`r' reg3b_`r' reg4_`r' reg4a_`r' reg4b_`r' using "centralitymeasure`r'.tex", label replace obslast star(* 0.10 ** 0.05 *** 0.01) se nomtitles mgroups("In-degree" "Out-degree" "Closeness" "Betweenness", pattern(1 0 0 1 0 0 1 0 0 1 0 0) prefix(\multicolumn{@span}{c}{) suffix(}) span erepeat(\cmidrule(lr){@span})) substitute(\_ _) ///
scalars("regressors Additional regressors$ \dagger $" "fixed Village/Enumerator FEs") b(2) sfmt(2) nonotes nogaps noomitted order(treatment treatmentmale blmale) drop(enumeratormale1 enumeratorfemale1 nightsaway readenglish hhsize distvillagemediankm farmer multiplespouses *.villagenumber)



local ++r

}



//////////////Robustness check///////////////

*Note: the code from this point takes approximately 10 hours to run on computer with specs: Intel Core i7 2.6 GHz, 32GB RAM

clear all
use "dyadicdatanetworkanalysis.dta", clear
sort villagenumber indvid indvidmch
drop if attrited==.
*total links: 25412
set seed 1

global bipred advice trust daily workedwith health closefriend
forvalues k=1/6{
gen random1=.
gen random2=.
gen rsample=.
local y=1
local l: word `k' of $bipred
forvalues x=.005(.005).20{
forvalues z=1/10{
replace random1=runiform()
sort random1
replace rsample=bipred_`l' if _n<=(`x'*25412)
putmata rsample, replace
replace random2=runiform()
sort random2

getmata rsample, replace

gen bipred_`l'robust_`y'_`z'=bipred_`l'
replace bipred_`l'robust_`y'_`z'=rsample if rsample!=.
replace bipred_`l'robust_`y'_`z'=. if indvid==indvidmch

}
local ++y
}
drop random1 random2 rsample
mata: mata clear
}
save "dyadicdatanetworkanalysisrobust.dta", replace
clear


local l=1
global outcomevars1 _in_degreestd _out_degreestd _closenessstd _betweenstd 
forvalues k=1/6{
matrix A=J(390,5,.)
local q=1

local l: word `k' of $bipred
forvalues x=1/39{
forvalues z=1/10{
foreach y in 1 2 3 5 7 8 10 11 12 14 16 17 18 19 21 24 26 27 28 30{
clear

global villagename Akibui Kongurakol Akiriamet Apokor Matisi Bokokholo Khasolo Abulei Emusioma Misikhu Malaha Sango Sirisia Naika Lunakwe Namisi Shipalo Shiraha Shirali Emuchelela Imanga Emokokho Milimani Kokare Kiriko Kocholia Lwanda EmukokhoA KapijanB Kakoit
matrix A[`q',1]=`x'

nwclear

use "dyadicdatanetworkanalysisrobust.dta", clear

keep if villagenumber==`y'
drop if attrited==.


nwset indvid indvidmch bipred_`l'robust_`x'_`z', name(network_indv_`x'_`z'_`y') edgelist keeporiginal 
rename _nodeoriginal indvid


local t: word `y' of $villagename

gen n=_n
egen maxn=max(n)
replace maxn=maxn-1

nwcloseness network_indv_`x'_`z'_`y', unconnected(max)

nwdegree network_indv_`x'_`z'_`y'


nwbetween network_indv_`x'_`z'_`y'
 

nwdegree network_indv_`x'_`z'_`y', valued
global outcomevars _in_degree _out_degree _closeness _between 

foreach var of varlist $outcomevars{
gen `var'std=.
sum `var' , detail 
replace `var'std=((`var'-`r(mean)')/`r(sd)') 
}

keep indvid _closenessstd _in_degreestd _out_degreestd _betweenstd 

save Stub`y', replace
}
nwclear
use Stub1, clear

foreach y in 2 3 5 7 8 10 11 12 14 16 17 18 19 21 24 26 27 28 30{
append using Stub`y'

}
save CentralityMeasures_`x'_`z', replace


use "alcdatacleanednets", clear
drop if indvid==.
drop if purecontrol==1
merge 1:1 indvid using CentralityMeasures_`x'_`z'
drop _merge
global villagecontrols blmale blage blalcoholnumber landareaacres assetindx nightsaway hhsize readenglish distvillagemediankm farmer yrseduc multiplespouses enumeratormale1 enumeratorfemale1 i.villagenumber
gen treatmentmale=treatment*blmale
gen controlmale=control*blmale
label var treatmentmale "Treatment assignment $ \times $ male"
label var distvillagemediankm "Distance to village center (KM)"

global villagenumber villagenum2 villagenum3 villagenum4 villagenum5 villagenum6 villagenum7 villagenum8 villagenum9 villagenum10 villagenum11 villagenum12 villagenum13 villagenum14 villagenum15 villagenum16 villagenum17 villagenum18 villagenum19 villagenum20 villagenum21 villagenum22 villagenum23 villagenum24 villagenum25 villagenum26 villagenum27 villagenum28 villagenum29 villagenum30
global outcomevars1 _in_degreestd _out_degreestd _closenessstd _betweenstd 

local u=2
foreach p of varlist $outcomevars1{

reg `p' treatment $villagecontrols, vce(bootstrap, reps(200) seed(2345) cluster(hhid))
gen t_stat`p'=_b[treatment]/_se[treatment]
matrix A[`q',`u']=t_stat`p'[1]

local ++u
}
local ++q
}
}
clear
svmat A
save "robust`l'", replace
egen maxA2=max(A2), by(A1)
egen minA2=min(A2), by(A1)
egen meanA2=mean(A2), by(A1)
egen maxA3=max(A3), by(A1)
egen minA3=min(A3), by(A1)
egen meanA3=mean(A3), by(A1)
egen maxA4=max(A4), by(A1)
egen minA4=min(A4), by(A1)
egen meanA4=mean(A4), by(A1)
egen maxA5=max(A5), by(A1)
egen minA5=min(A5), by(A1)
egen meanA5=mean(A5), by(A1)
collapse maxA2 minA2 meanA2 maxA3 minA3 meanA3 maxA4 minA4 meanA4 maxA5 minA5 meanA5, by(A1)
drop if maxA2==.
twoway rbar maxA2 minA2 A1, ytitle("T-statistic on Treatment") xlabel(0 "0" 9 "5" 19 "10" 29 "15" 39 "20") xtitle("Percent of randomly replaced links") legend(label(1 "Estimate Range") label(2 "Mean value")) || scatter meanA2 A1, msymbol(circle) mcolor(black) msize(tiny)
graph export "robust`l'_indegreegraph.png", replace
twoway rbar maxA3 minA3 A1, ytitle("T-statistic on Treatment") xlabel(0 "0" 9 "5" 19 "10" 29 "15" 39 "20") xtitle("Percent of randomly replaced links") legend(label(1 "Estimate Range") label(2 "Mean value")) || scatter meanA3 A1, msymbol(circle) mcolor(black) msize(tiny)
graph export "robust`l'_outdegreegraph.png", replace
twoway rbar maxA4 minA4 A1, ytitle("T-statistic on Treatment") xlabel(0 "0" 9 "5" 19 "10" 29 "15" 39 "20") xtitle("Percent of randomly replaced links") legend(label(1 "Estimate Range") label(2 "Mean value")) || scatter meanA4 A1, msymbol(circle) mcolor(black) msize(tiny)
graph export "robust`l'_closegraph.png", replace
twoway rbar maxA5 minA5 A1, ytitle("T-statistic on Treatment") xlabel(0 "0" 9 "5" 19 "10" 29 "15" 39 "20") xtitle("Percent of randomly replaced links") legend(label(1 "Estimate Range") label(2 "Mean value")) || scatter meanA5 A1, msymbol(circle) mcolor(black) msize(tiny)
graph export "robust`l'_betweengraph.png", replace

}
